Для выполнения этого задания вам понадобятся данные о кредитных историях клиентов одного из банков. Поля в предоставляемых данных имеют следующий смысл:
In [1]:
from __future__ import division
import numpy as np
import pandas as pd
from scipy import stats
from statsmodels.stats import weightstats as wsm
from statsmodels.stats.proportion import proportion_confint
%matplotlib inline
import matplotlib.pyplot as plt
import seaborn as sns
from IPython.core.interactiveshell import InteractiveShell
InteractiveShell.ast_node_interactivity = "all"
In [2]:
#reading the data
credit_story = pd.read_csv('credit_card_default_analysis.csv')
credit_story.info()
In [3]:
pd.set_option('display.max_columns', None)
credit_story.head()
Out[3]:
LIMIT_BAL: размер кредитного лимита (в том числе и на семью клиента)
In [4]:
#histogram of LIMIT_BAL distribution
_ = plt.figure(1, figsize=(15,7))
_ = plt.subplot(121)
_ = plt.title('Histogram of LIMIT_BALL')
_ = plt.hist(credit_story['LIMIT_BAL'], edgecolor='k')
_ = plt.subplot(122)
_ = plt.title('Histogram of LIMIT_BALL for default = 0 / 1')
_ = plt.hist(credit_story.LIMIT_BAL.loc[credit_story.default == 0], edgecolor='k', label='Default = 0')
_ = plt.hist(credit_story.LIMIT_BAL.loc[credit_story.default == 1], edgecolor='k', label='Default = 1')
_ = plt.legend()
Доля небольших невозвратных кредитов заметно выше. Большие суммы люди стремятся возвращать.
In [5]:
#samples of LIMIT_BAL for default = 0/1
lim_bal_0 = credit_story.LIMIT_BAL.loc[credit_story.default == 0]
lim_bal_1 = credit_story.LIMIT_BAL.loc[credit_story.default == 1]
print('Median LIMIT_BAL for default = 0: %.0f' % lim_bal_0.median())
print('Median LIMIT_BAL for default = 1: %.0f' % lim_bal_1.median())
In [6]:
def zconfint_binom(n, p, alpha=0.05):
q = 1 - p
m = n * p
var = n * p * q
z_stat = stats.norm.ppf(1 - alpha / 2)
l_bound = int(round(m - z_stat * np.sqrt(var)))
u_bound = int(round(m + z_stat * np.sqrt(var)))
return (l_bound, u_bound)
In [7]:
#confidence intervals estimation
lim_bal_0_s = np.sort(lim_bal_0.values)
lim_bal_1_s = np.sort(lim_bal_1.values)
l_ind_0, u_ind_0 = zconfint_binom(len(lim_bal_0_s), 0.5)
l_ind_1, u_ind_1 = zconfint_binom(len(lim_bal_1_s), 0.5)
print('Conf. int. median LIMIT_BAL for default = 0: [%d, %d]' % (lim_bal_0_s[l_ind_0], lim_bal_0_s[u_ind_0]))
print('Conf. int. median LIMIT_BAL for default = 1: [%d, %d]' % (lim_bal_1_s[l_ind_1], lim_bal_1_s[u_ind_1]))
In [8]:
def get_bootstrap_samples(data, n_samples):
indices = np.random.randint(0, len(data), (n_samples, len(data)))
samples = data[indices]
return samples
In [9]:
def stat_intervals(stat, alpha):
boundaries = np.percentile(stat, [100 * alpha / 2., 100 * (1 - alpha / 2.)])
return boundaries
In [10]:
#confidence intervals estimation
np.random.seed(0)
lim_bal_0_bs_med = map(np.median, get_bootstrap_samples(lim_bal_0.values, 100))
lim_bal_1_bs_med = map(np.median, get_bootstrap_samples(lim_bal_1.values, 100))
bnd_0 = stat_intervals(lim_bal_0_bs_med, 0.05)
bnd_1 = stat_intervals(lim_bal_1_bs_med, 0.05)
print('Conf. int. median LIMIT_BAL for default = 0: [%d, %d]' % (bnd_0[0], bnd_0[1]))
print('Conf. int. median LIMIT_BAL for default = 1: [%d, %d]' % (bnd_1[0], bnd_1[1]))
Как видно из интервальной оценки, медианы не совпадают. Заёмщики, которые не возвращают кредит, обычно берут меньшую сумму.
Для проверки гипотезы подойдёт перестановочный критерий для независимых выборок.
In [11]:
def permutation_t_stat_ind(sample1, sample2):
return np.mean(sample1) - np.mean(sample2)
In [12]:
def get_random_combinations(n1, n2, max_combinations):
index = range(n1 + n2)
indices = set([tuple(index)])
for i in range(max_combinations - 1):
np.random.shuffle(index)
indices.add(tuple(index))
return [(index[:n1], index[n1:]) for index in indices]
In [13]:
def permutation_zero_dist_ind(sample1, sample2, max_combinations=None):
joined_sample = np.hstack((sample1, sample2))
n1 = len(sample1)
n = len(joined_sample)
if max_combinations:
indices = get_random_combinations(n1, len(sample2), max_combinations)
else:
indices = [(list(index), filter(lambda i: i not in index, range(n))) \
for index in itertools.combinations(range(n), n1)]
distr = [joined_sample[list(i[0])].mean() - joined_sample[list(i[1])].mean() \
for i in indices]
return distr
In [14]:
def permutation_test(sample1, sample2, max_permutations=None, alternative='two-sided'):
if alternative not in ('two-sided', 'less', 'greater'):
raise ValueError("alternative not recognized\n"
"should be 'two-sided', 'less' or 'greater'")
t_stat = permutation_t_stat_ind(sample1, sample2)
zero_distr = permutation_zero_dist_ind(sample1, sample2, max_permutations)
if alternative == 'two-sided':
return sum([1. if abs(x) >= abs(t_stat) else 0. for x in zero_distr]) / len(zero_distr)
if alternative == 'less':
return sum([1. if x <= t_stat else 0. for x in zero_distr]) / len(zero_distr)
if alternative == 'greater':
return sum([1. if x >= t_stat else 0. for x in zero_distr]) / len(zero_distr)
In [15]:
print('p-value: %f' % permutation_test(lim_bal_0, lim_bal_1, max_permutations=1000))
P-value получается довольно маленьким и нулевая гипотеза отвергается на уровне значимости 0.05.
Результат является практически значимым, потому что разница в значениях медианы LIMIT_BAL значимо большая.
SEX: пол клиента (1 = мужской, 2 = женский )
In [16]:
#histogram of SEX distribution
_ = plt.figure(1, figsize=(15,7))
_ = plt.subplot(121)
_ = plt.title('Histogram of SEX')
_ = plt.hist(credit_story['SEX'], edgecolor='k')
_ = plt.subplot(122)
_ = plt.title('Histogram of SEX for default = 0 / 1')
_ = plt.hist(credit_story.SEX.loc[credit_story.default == 0], edgecolor='k', label='Default = 0')
_ = plt.hist(credit_story.SEX.loc[credit_story.default == 1], edgecolor='k', label='Default = 1')
_ = plt.legend()
Женщины более склонны к невозврату кредита, чем мужчины.
In [17]:
#samples of SEX for default = 0/1
sex_0 = credit_story.SEX.loc[credit_story.default == 0].values
sex_1 = credit_story.SEX.loc[credit_story.default == 1].values
#proportions of men to women in samples
m_to_w_0 = np.where(sex_0 == 1)[0].shape[0] / sex_0.shape[0]
m_to_w_1 = np.where(sex_1 == 1)[0].shape[0] / sex_1.shape[0]
print('Prop. of men for default = 0: %.4f' % m_to_w_0)
print('Prop. of men for default = 1: %.4f' % m_to_w_1)
Воспользуемя доверительным интервалом на основе нормального распределения.
In [18]:
#confident intervals for proportions
conf_int_0 = proportion_confint(np.where(sex_0 == 1)[0].shape[0], sex_0.shape[0])
conf_int_1 = proportion_confint(np.where(sex_1 == 1)[0].shape[0], sex_1.shape[0])
print('Conf. int. for prop. of men or default = 0: [%.4f, %.4f]' % conf_int_0)
print('Conf. int. for prop. of men or default = 1: [%.4f, %.4f]' % conf_int_1)
Как видно из интервальных оценок, сами интервалы не пересекаются, следовательно гендерный состав различается.
Построим доверительный интервал для разности двух долей независимых выборок.
In [19]:
def proportions_confint_diff_ind(sample1, sample2, alpha = 0.05):
z = stats.norm.ppf(1 - alpha / 2.)
p1 = float(sum(sample1)) / len(sample1)
p2 = float(sum(sample2)) / len(sample2)
left_boundary = (p1 - p2) - z * np.sqrt(p1 * (1 - p1) / len(sample1) + p2 * (1 - p2) / len(sample2))
right_boundary = (p1 - p2) + z * np.sqrt(p1 * (1 - p1) / len(sample1) + p2 * (1 - p2) / len(sample2))
return (left_boundary, right_boundary)
In [20]:
print('Conf. int. for difference: [%f, %f]' % proportions_confint_diff_ind(np.abs(sex_0-2), np.abs(sex_1-2)))
In [21]:
def proportions_diff_z_stat_ind(sample1, sample2):
n1 = len(sample1)
n2 = len(sample2)
p1 = float(sum(sample1)) / n1
p2 = float(sum(sample2)) / n2
P = float(p1*n1 + p2*n2) / (n1 + n2)
return (p1 - p2) / np.sqrt(P * (1 - P) * (1. / n1 + 1. / n2))
In [22]:
def proportions_diff_z_test(z_stat, alternative='two-sided'):
if alternative not in ('two-sided', 'less', 'greater'):
raise ValueError("alternative not recognized\n"
"should be 'two-sided', 'less' or 'greater'")
if alternative == 'two-sided':
return 2 * (1 - stats.norm.cdf(np.abs(z_stat)))
if alternative == 'less':
return stats.norm.cdf(z_stat)
if alternative == 'greater':
return 1 - stats.norm.cdf(z_stat)
In [23]:
print('p-value: %.15f' % proportions_diff_z_test(proportions_diff_z_stat_ind(np.abs(sex_0-2), np.abs(sex_1-2))))
Гипотеза о равенстве гендерного распределния отвергается. Данный результат явлется и практически значимым, т.к. 6% существенная величина для банков.
EDUCATION: образование (0 = доктор, 1 = магистр; 2 = бакалавр; 3 = выпускник школы; 4 = начальное образование; 5 = прочее; 6 = нет данных)
In [24]:
#histogram of EDUCATION distribution
_ = plt.figure(1, figsize=(15,7))
_ = plt.subplot(121)
_ = plt.title('Histogram of EDUCATION')
_ = plt.hist(credit_story['EDUCATION'], edgecolor='k')
_ = plt.subplot(122)
_ = plt.title('Histogram of EDUCATION for default = 0 / 1')
_ = plt.hist(credit_story.EDUCATION.loc[credit_story.default == 0], edgecolor='k', label='Default = 0')
_ = plt.hist(credit_story.EDUCATION.loc[credit_story.default == 1], edgecolor='k', label='Default = 1')
_ = plt.legend()
В основном берут кредиты люди с BS и MS. Принципиальная разница между распределниями тех, кто вернул, и тех, кто нет, не наблюдается.
Построим график долей возврата кредитов.
In [25]:
#samples of EDUCATION for default = 0/1
edu = credit_story.EDUCATION
edu_0 = credit_story.EDUCATION.loc[credit_story.default == 0]
edu_1 = credit_story.EDUCATION.loc[credit_story.default == 1]
#proportions of credit returns
edu_prop = np.empty( (len(np.unique(edu_0)), 2) )
for i, ed_val in enumerate(np.unique(edu)):
edu_prop[i, 0] = ed_val
edu_prop[i, 1] = edu_0[edu_0 == ed_val].shape[0] / (edu_0[edu_0 == ed_val].shape[0] + edu_1[edu_1 == ed_val].shape[0])
for ed_val, prop in edu_prop:
print('Prop. of def.=0 for education=%d: %.2f' % (ed_val, prop))
#histogram of EDUCATION proportions distribution
_ = plt.figure(1, figsize=(7,7))
_ = plt.title('Proportion histogram of EDUCATION')
_ = plt.plot(edu_prop[:,0], edu_prop[:,1])
_ = plt.axis([0, 6, 0.5, 1])
Исходя из графика видно, что распределние не является равномерным, следовательно, образование влияет на то, вернёт ли человек долг. Реже всего возвращают долги категории 2 и 3 (бакалавры и выпускники школ). Возвращают практически все долги только люди с докторской степенью.
Посчитаем среднее значение доли возврата кредитов среди различных уровней образования. Будем использовать эту оценку для построения ожидаемого равномерного распределения.
In [26]:
prop_mean = edu_prop[:,1].mean()
print('Mean proportion value: %.2f' % prop_mean)
In [27]:
#sample sizes for different education levels
edu_num_0 = np.array([edu_0[edu_0 == ed_val].shape[0] for ed_val in range(7)])
edu_num_1 = np.array([edu_1[edu_1 == ed_val].shape[0] for ed_val in range(7)])
edu_num = edu_num_0 + edu_num_1
for i in range(7):
print('Edu lvl %d: %d / %d | %d' % (i, edu_num_0[i], edu_num_1[i], edu_num[i]))
print('Legend: deafault = 0 / default = 1 | total')
Составим таблицу сопряжённости "образование" на "возврат долга", где значением ячейки была бы разность между наблюдаемым и ожидаемым количеством человек.
In [28]:
#expected frequencies
exp_freq = np.array([prop_mean * edu_num[i] for i in range(7)])
for i in range(7):
print('Edu lvl %d: %d / %d | %d' % (i, exp_freq[i], edu_num_0[i], edu_num_0[i] - exp_freq[i]))
print('Legend: exp freq / obs freq = 1 | delta')
In [29]:
#chi-square criterion
stats.chisquare(edu_num_0, exp_freq, ddof = 1)
Out[29]:
Гипотеза о том, что величина имеет равномерное распределение отвергается с очень высокой долей вероятности.
Результат имеет практическую значимость для выдачи кредита, доли возврата значимо отличаются для разных уровней образования.
MARRIAGE: (0 = отказываюсь отвечать; 1 = замужем/женат; 2 = холост; 3 = нет данных).
In [30]:
#histogram of MARRIAGE distribution
_ = plt.figure(1, figsize=(15,7))
_ = plt.subplot(121)
_ = plt.title('Histogram of MARRIAGE')
_ = plt.hist(credit_story['MARRIAGE'], edgecolor='k')
_ = plt.subplot(122)
_ = plt.title('Histogram of MARRIAGE for default = 0 / 1')
_ = plt.hist(credit_story.MARRIAGE.loc[credit_story.default == 0], edgecolor='k', label='Default = 0')
_ = plt.hist(credit_story.MARRIAGE.loc[credit_story.default == 1], edgecolor='k', label='Default = 1')
_ = plt.legend()
Доли невозвратов среди женатых и холостых людей различается примерно на 10%.
In [31]:
#samples of MARRIAGE for default = 0/1
mar = credit_story.MARRIAGE
mar_0 = credit_story.MARRIAGE.loc[credit_story.default == 0]
mar_1 = credit_story.MARRIAGE.loc[credit_story.default == 1]
print('Value counts for marriage:')
mar.value_counts()
Out[31]:
Составим таблицу сопряжённости для признаков default и MARRIAGE.
In [32]:
conf_table = np.empty( (len(mar.value_counts()), 2) )
for i in range(len(mar.value_counts())):
conf_table[i, 0] = len(mar_0.loc[mar_0 == i])
conf_table[i, 1] = len(mar_1.loc[mar_1 == i])
conf_table
Out[32]:
Вычислим коэффициент V Крамера.
In [33]:
def cramers_stat(confusion_matrix):
chi2 = stats.chi2_contingency(confusion_matrix)[0]
n = confusion_matrix.sum()
return np.sqrt(chi2 / (n*(min(confusion_matrix.shape)-1)))
In [34]:
print('V Cramer stat value: %.4f' % cramers_stat(conf_table))
Коэффициент Крамера принимает значение, близкое к нулю, когда взаимосвязь отсутствует.
Практическая значимость в данном результате также присутствует, т.к. исходя из этих данных нет разницы в том, с каким статусом MARRIAGE приходит человек за кредитом.
In [35]:
#histogram of AGE distribution
_ = plt.figure(1, figsize=(15,7))
_ = plt.subplot(121)
_ = plt.title('Histogram of AGE')
_ = plt.hist(credit_story['AGE'], edgecolor='k')
_ = plt.subplot(122)
_ = plt.title('Histogram of AGE for default = 0 / 1')
_ = plt.hist(credit_story.AGE.loc[credit_story.default == 0], edgecolor='k', label='Default = 0')
_ = plt.hist(credit_story.AGE.loc[credit_story.default == 1], edgecolor='k', label='Default = 1')
_ = plt.legend()
Молодые люди склонные не возвращать кредиты. Доля тех, кто взял кредит в 20 лет и вернул его, меньше доли оных, которые не вернули кредит.
In [36]:
#samples of AGE for default = 0/1
age_0 = credit_story.AGE.loc[credit_story.default == 0]
age_1 = credit_story.AGE.loc[credit_story.default == 1]
print('Median AGE for default = 0: %.0f' % age_0.median())
print('Median AGE for default = 1: %.0f' % age_1.median())
In [37]:
#confidence intervals estimation
np.random.seed(0)
age_0_bs_med = map(np.median, get_bootstrap_samples(age_0.values, 100))
age_1_bs_med = map(np.median, get_bootstrap_samples(age_1.values, 100))
bnd_0 = stat_intervals(age_0_bs_med, 0.05)
bnd_1 = stat_intervals(age_1_bs_med, 0.05)
print('Conf. int. median AGE for default = 0: [%f, %f]' % (bnd_0[0], bnd_0[1]))
print('Conf. int. median AGE for default = 1: [%f, %f]' % (bnd_1[0], bnd_1[1]))
Как видно из интервальной оценки, интервальные оценки для медианы пересекаются.
In [38]:
print('Difference between medians: %f' % (np.median(age_0) - np.median(age_1)))
In [39]:
delta_median_scores = map(lambda x: x[1] - x[0], zip(age_0_bs_med, age_1_bs_med))
In [40]:
stat_int = stat_intervals(delta_median_scores, 0.05)
print('95%% conf. int. for the difference between medians: [%f, %f]' % (stat_int[0], stat_int[1]))
Поскольку доверительный интервал содержит ноль, значения медиан не различаются.
Для проверки гипотезы подойдёт перестановочный критерий для независимых выборок.
In [41]:
print('p-value: %f' % permutation_test(age_0, age_1, max_permutations=1000))
P-value получается довольно маленьким и нулевая гипотеза отвергается на уровне значимости 0.05.
Результат является практически значимым, потому что разница в распределениях видна невооружённым взглядом и эти данные могут учитываться при выдаче кредита.